Glioblastoma pseudoprogression and true progression reveal spatially variable transcriptional differences

Post-resection radiologic monitoring to identify areas of new or progressive enhancement concerning for cancer recurrence is critical during patients with glioblastoma follow-up. However, treatment-related pseudoprogression presents with similar imaging features but requires different clinical management. While pathologic diagnosis is the gold standard to differentiate true progression and pseudoprogression, the lack of objective clinical standards and admixed histologic presentation creates the needs to (1) validate the accuracy of current approaches and (2) characterize differences between these entities to objectively differentiate true disease. We demonstrated using an online RNAseq repository of recurrent glioblastoma samples that cancer-immune cell activity levels correlate with heterogenous clinical outcomes in patients. Furthermore, nCounter RNA expression analysis of 48 clinical samples taken from second neurosurgical resection supports that pseudoprogression gene expression pathways are dominated with immune activation, whereas progression is predominated with cell cycle activity. Automated image processing and spatial expression analysis however highlight a failure to apply these broad expressional differences in a subset of cases with clinically challenging admixed histology. Encouragingly, applying unsupervised clustering approaches over our segmented histologic images provides novel understanding of morphologically derived differences between progression and pseudoprogression. Spatially derived data further highlighted polarization of myeloid populations that may underscore the tumorgenicity of novel lesions. These findings not only help provide further clarity of potential targets for pathologists to better assist stratification of progression and pseudoprogression, but also highlight the evolution of tumor-immune microenvironment changes which promote tumor recurrence. Supplementary Information The online version contains supplementary material available at 10.1186/s40478-023-01587-w.

Table S2: Admixed sample lesion description and pathology.Relevant clinical information about the 8 samples used in both our image analysis and DSP studies are shown in the table.All samples were confirmed by pathology to have a mixed presentation as listed in the pathology diagnostic report.Relevant patient protected health information has been removed.Time to enhancement was listed based on nearest months from ChemoRT completion to novel enhancement detection on radiologic imaging.Slide tissue area refers to the area of tissue present on the slide used for DSP analysis.Lesion size was listed by the 2 largest axes measured on clinical report.Ki67% was measured by pathologic interpretation taken from the pathology report.Side, lobe, and diagnosis were all collected from the neuro-oncology report.
Table S3: Mixed effect modeling statistics of image analysis groupings to segmentation measures.Mixed effect ANOVA scores were calculated and reported in the table for novel enhancement status and combined status plus histology.Relevant values for ANOVA score are listed in numerator degrees of freedom (NumDF) versus denominator degrees of freedom (DenDF), F-value, and P-value.Effect size for each grouping is reported by Eta-squared whereby larger scores represent greater effect on the measured outcome.Significance with mixed ANOVA was reported in the main figures as a red asterisk.Calculation of scores was performed in R using linear mixed method calculation to control on random effect due to multiple sample testing.
Table S4: Mixed effect modeling statistics of DSP analysis groupings to WGCNA and immune deconvolution measures.Mixed effect ANOVA scores were calculated and reported in the table for novel enhancement status and combined status plus histology.Relevant values for ANOVA score are listed in numerator degrees of freedom (NumDF) versus denominator degrees of freedom (DenDF), F-value, and P-value.Effect size for each grouping is reported by Etasquared whereby larger scores represent greater effect on the measured outcome.Significance with mixed ANOVA was reported in the main figures as a red asterisk.Calculation of scores was performed in R using linear mixed method calculation to control on random effect due to multiple sample testing.

Figure S1
: Tissue processing framework for molecular studies.FFPE tissues were processed and evaluated using bulk nCounter analysis in the upper path.Tissue was extracted for RNA, measured for RNA integrity and concentration, and processed on a nCounter Pancancer 360 panel.Tissue was sectioned and imaged for GeoMx DSP in the lower path.Tissue was sectioned at 5 microns and stained with conjugated antibodies and probes.Regions of interest were drawn on GeoMx to cleave oligo-tags.Final ground truth designation for binning of cases is based on clinico-pathologic diagnosis rendered by neuro-oncology.Full details of processing can be found in the methods, DEGs between PD and psPD GB novel enhancement events shown in nCounter studies.Left-side hierarchical clustering had grouped similar expressing genes with names listed on the right.Topside hierarchical clustering had grouped similar samples with PD in teal or psPD in red.Specific sample label is found on the bottom of the heatmap.kMeans clustering in main figure 2 is independent of the clustering in the heatmap; however, same genes were used in the analysis.Highlighting in either case that not all samples perfectly stratify, Figure S3: Image processing schematic of IHC stains.(A) Segmentation and analysis scheme for H&E images [n=260 for 8 cases].Images were deconvoluted and segmented to calculate morphology features.Prior to performing either nearest neighbor spatial distribution calculation or morphology clustering, segmentations were filtered based on size to remove punctate objects.(B) Segmentation and analysis scheme for Olig2, Ki67, and p53 images [n=260 per stain for 8 cases].Images were deconvoluted and segmented to calculate morphology features.However, a random forest machine classifier was trained to identify segments as true signal or not based on positive and negative controls.This was due to the variable staining intensity across images.After filtration, pixel stain ratio calculation and morphology clustering were done.(C) CD163 scheme for (1) segmentation and analysis of image tiles [n=5,809 for 8 cases] and (2) CNN model training [n=6,524 images].Development of CNN models was done with an external dataset of CD163 images that were manually annotated by a trained observer.Representation of brown background stain was done by computationally noising the image with brown pixels.After training, images taken from the 8 samples were passed through both models to generate prediction maps of cell processes and bodies.This was done due to the different optimizations needed to accurately identify processes and bodies.Images were then remerged to calculate morphology features and filtered for downstream processing as described in B. (D) Intersection over union accuracy score of CNN models after mask merge against ground truth test dataset [n=1,288].IoU score was based on accurate pixel categorization of both negative and positive pixels as opposed to an accuracy score.This is based on the fact that accuracy would measure the total number of pixels correctly guessed as positive without integrating negative pixels-causing a bias as most images do not have many positive pixels.Outliers with low score are shown at the bottom of the boxplot-showing discordance of ground truth image and CNN prediction.(E) Representation of outlier segmentations in IoU score showing inaccuracy due to lack of segmentations in ground truth.Ground truths in external dataset were created by a trained observer manually.In consequence, missed segmentations are possible in the ground truth images that the CNN overcame.

Figure S4:
Cluster statistics of CGGA clustering analysis.Clustering statistics of (A) average within cluster distance and (B) average between cluster distance of points in actual dataset (green) and randomly scrambled dataset (grey).Randomly scrambled dataset was generated by scrambling values taken from the actual data-effectively eliminating potential relationships in the dataset.Average within cluster distance represents the mean distance a point is from all other points in a cluster, while average between cluster distance represent mean distance of a point from all points in another cluster.Statistically valid clustering would be shown as significant deviation of clustered point scores against a scrambled dataset where no real relationship is present.T-test: **** p ≤ 0.0001 Figure S5: WGCNA module and immune proliferation gene overlap.Overlap of genes related to lymphocyte proliferation (GO:0046651) [lower left] or macrophage proliferation (GO:0061517) [lower right] for CGGA WGCNA color modules related to (A) cell cycle activity and (B) immune processes.Gene names were extracted from all found groups and compared using a Venn diagram to identify whether overlapped use of terms was present.In specific, query was done to see whether our color modules contained genes affiliated to immune proliferation, Figure S6: Cluster statistics of nCounter clustering analysis.Clustering statistics of (A) average within cluster distance and (B) average between cluster distance of points in actual dataset (green) and randomly scrambled dataset (grey).Randomly scrambled dataset was generated by scrambling values taken from the actual data-effectively eliminating potential relationships in the dataset.Average within cluster distance represents the mean distance a point is from all other points in a cluster, while average between cluster distance represent mean distance of a point from all points in another cluster.Statistically valid clustering would be shown as significant deviation of clustered point scores against a scrambled dataset where no real relationship is present.T-test: **** p ≤ 0.0001 Figure S7: Representative staining of IHC markers in sub-stratified samples.Representative imaging of IHC stains taken from our 8 admixed samples to highlight variation in PD and psPD events in spite of similar overall histology categorization (control, hypercellular, and inflammatory).Samples are captured at 40x resolution from Phillips Image Management System.It can be appreciated that overall histology between events is similar, but subtle differences in staining intensity and morphology are present.
Figure S8: Volcano plot analysis of differentially expressed genes in "control" regions in PD and psPD event from GeoMx.Differential analysis of PD (right) and psPD (left) events with respect to "control" regions are highlighted in the volcano plot.Labeled genes represent those that show significant log fold change towards a direction of disease."Control" histology is represented as overall normal brain with typical neurons and oligodendrocytes present.

Figure S2 :
Figure S2: Top 150 DEGs from OSU PD vs psPD Cases.Expression heatmap of top 150DEGs between PD and psPD GB novel enhancement events shown in nCounter studies.Left-side hierarchical clustering had grouped similar expressing genes with names listed on the right.Topside hierarchical clustering had grouped similar samples with PD in teal or psPD in red.Specific sample label is found on the bottom of the heatmap.kMeans clustering in main figure2is independent of the clustering in the heatmap; however, same genes were used in the analysis.Highlighting in either case that not all samples perfectly stratify,